Mapping the Two-Component Atomic Fermi Gas to the Nuclear Shell-Model 
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The physics of two-component fermi systems is now frequently addressed in laboratories. Usually 
this is done for large samples of tens to hundreds of thousands of particles. However, there are 
now groups working towards few-body systems (1-100 particles) in very tight traps where the shell 
structure of the external potential becomes important. This is analogous to the situation encountered 
in nuclear structure physics. Here the state-of-the-art method is large-scale diagonalization of the 
interaction hamiltonian in an appropriate model space. We map the atomic physics problem onto the 
nuclear shell-model and show that the two-component atomic fermi gas with zero-range interaction 
does not suffer from the famous fermionic sign problem, making Quantum Monte Carlo methods 
feasible. 



I. INTRODUCTION 



The physics of ultracold gases has seen rapid developments over the past two decades [T|, 0, Q ■ One of the next 
goals is to trap few fermions to approach the regime where the trapping shell structure becomes prominent Q. For 
small number of particles one expects large quantum fluctuations that invalidate the use of mean-field theories such 
as the BCS approach. Theoretical investigations beyond mean-field are therefore needed and on-going [!, H[. In this 
paper we demonstrate how methods to solve the many-body problem in nuclei can be applied by mapping the atomic 
two-component fermion system in a tight trap into that of a single species of nucleon in an external field. 

One particular successful approach to the nuclear shell- model is the so-called Shell-model Monte Carlo method 0] • 
lO \ This method overcomes the large-scale matrix diagonalization by transforming two-body terms into one-body terms 
04 ■ and subsequently solves a multi-dimensional path integral using Monte Carlo integration. It has had many successes 
in nuclear application and it is hoped that one might adapt this method to also address small atomic fermi systems. 
However, there is an infamous sign problem that plagues Monte Carlo methods when dealing with fermions [8J, and 
s ! ' it is not a priory obvious that the atomic two-component Fermi system can be handled at all. 

In this article we discuss this mapping in some detail and give the essentials of the shell-model Monte Carlo method. 
We then proceed to prove the absence of a sign problem for the case of a zero-range interaction. We also briefly discuss 
of some similarities and differences of pairing studies in atomic and nuclear physics. For the benefit of both nuclear 
and atomic physics communities, the explicit evaluation of the zero-range interaction relevant for the atomic physics 
system in terms of the matrix elements used in the nuclear shell model is given in considerable detail in an appendix. 



II. MAPPING THE FERMI GAS TO THE NUCLEAR SHELL MODEL 



The two-component ultracold fermionic atomic gas consists of neutral atoms, usually alkali ones, that occupy two 
different internal states. The actual internal states are hyperfine states of different projection that can be split by a 
magnetic field Q. The hyperfine splitting is larger than any other energy scale in the problem so that no internal 
process in the gas can transfer atoms between the hyperfine levels, thus one may think of these as frozen degrees of 
freedom. Also, since the systems are usually dilute, the range of the atom-atom interactions is very short compared to 
the typical interparticle distance. Therefore the simple zero-range potential is a popular and highly successful model. 

The ./V-body fermi system in an isotropic harmonic trap with a zero-range interaction of strength Vq has the 
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hamiltonian 

^ = E|J+E + E w< - d) 

» * [*i] 

where [ij] denotes the sum over all pairs of particles. In terms of two-body matrix elements we have 

{ip a (n ) (fa ) I V mt \^c (n (fa)) , (2) 

where Vi„t is any two-body interaction. The two-body wavefunction V'^f'OV'fc^) and ^0(^1)^(^2) must be antisym- 
metric under exchange of the coordinates. 

In the tight harmonic trap, the quantum numbers are n,l,mi. The two internal hyperfine states can now be 
mapped onto the spin of a single species of nucleon m s = ±1/2. Thus, a single-particle state is uniquely given 
by a = (n a l a mi a m Sa ). Any two-body state constructed from these states will have external and internal parts 
that combine to determine the overall symmetry. Since we use a zero-range interaction, the spatial part must be 
symmetric under particle exchange to give a non-zero contribution. Therefore the internal (pseudospin) state must 
be antisymmetric, which in the spin 1/2 language is the spin-singlet state. This completes our mapping of the two- 
component fermi system in a trap onto the nuclear shell model with one species of nucleon (proton or neutron). The 
nuclear mean-field is simply replaced by the harmonic oscillator and the internal spin states of the nucleon are now 
hyperfine states for the atoms. 



III. AUXILIARY-FIELD MONTE CARLO METHOD 

The use of Quantum Monte Carlo for fermions is plagued by the well-known sign problem (see f.x. the review 0]). 
This basically arises because there is no way to choose a positive definite measure for the Monte Carlo sampling when 
taking into account the exchange symmetry of fermions. In fact it was shown in [3| that the problem is computationally 
NP hard. Various ways have been invented to deal with the problem by appropriately restricting the Hilbert space or 
interacting hamiltonian of the system under study. One possibility is to integrate out the fermions and get an effective 
bosonic theory [lol [ll[ , another is to use variational wavefunctions with a fixed node structure (T3 , [Tsl , a third way 
is the so-called determinant diagrammatic Monte Carlo which is also principally sign problem free [14j . Common to 
all the methods above is that they are mostly successful for homogeneous fermi systems. Here we are concerned with 
a very different limit of small systems in tight traps, thus the sensible quantum numbers are quite different. This 
corresponds to the situation in nuclear physics. 

Here we review the Auxiliary-field Monte Carlo method for the fermionic many-body problem. In the context of 
the nuclear shell- model this approach is known as the Shell-Model Monte Carlo (SMMC) method [7(. This method is 
based on a linearization of the two-body part of the Hamiltonian using the Hubbard-Stratonovich transformation. Here 
we adopt a formulation of this transformation starting from a general Hamiltonian which is explicitly time-reversal 
invariant: 

H = ^(e a O a + e*0 a ) + i^7 a {O a ,O a } ) (3) 

a a 

where O a are one-body operators in a convenient basis and the V a are real numbers. The bars denote time-reversed 
operators. The SMMC approach uses the Hubbard-Stratonovich (HS) transformation to linearize the many-body 
evolution operator e~@ H , where /3 _1 may be interpreted as the temperature in the (grand)canonical ensemble. We 

first divide (3 into Nt time slices so that we can express individual terms at different time slices in e~P H = [e~ A @ H ] Nt 
as 

a 

where we used 2 {O a ,O a } = (O a + O a ) 2 — (O a — O a ) 2 . Quadratic interaction terms can be effectively linearized 
through the Gaussian integral identity 



da^da a e 



-A/3^[(<^) 2 + (0 2 



47T 

e ~A/3^[s a cr^(O a +O a )+is a <T I a (O a ~O a )] 



(4) 
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where the integration variables erf and a T a are the real auxiliary fields that give the method its name. The sign 
factors are s a = ±1 for V a < and s a — ±z for V a > 0. Introducing complex fields for each time slice cr a ( T n) = 
a a( T n) + io^(r n ), we arrive at the Hubbard-Stratonovich representation of the many-body evolution operator 



e -PH = 



V[a]G{a)U a (/3,0). (5) 



Above, 



jyy\ _ TT dgq(7-n)^(r n ) A^pVg 
■LA 2i 47r 

is the measure of the integral. G(a) is a Gaussian weight 

Many-body propagator is now effectively reduced to a superposition of one-body propagators 

U„(p,0) =e -A^«(r Wt ) > __ e -A^ CT (r 1 ) ) (g) 

where the linearized Hamiltonian as a function of the time-dependent auxiliary fields is given by 

h °( T ) = E ( e « + ^ s «^«(t)J O q + (e* a + \s a V a a* a {T)^j 6 a . (9) 

In the SMMC, expectation value of an observable ft is calculated by expressing both the numerator and the de- 
nominator of (ft) = Trjv[f2e -,3ff ]/Trjve~^ ff (where Tr^y denote canonical traces for N-particle system) in the HS 
representation. In order to perform a Monte Carlo integration, a positive definite weight function is defined as 
W(a) = G(<r)\TipfU a (P, 0)|. Thus, one can express the thermal expectation values by 

JV[a}W(aWa) (10) 

where $(cr) = Tr N U a {P, 0)/|Tr JV [/ (7 (/3, 0)| is the "sign" and = Tr N \nU a (/3,0)\/Tr N U a (/3,0). The observable 

(Q) is then computed in a Monte Carlo integration by selecting an ensemble of auxiliary fields (<ti, . . . , ajq) sampled 
according to the distribution function W(a), i.e., 

(ft) =w ^n%)("). (n) 

Success of the outlined method hinges on the sign Q(a) of the weight function W(a). Unfortunately, in the most general 
case, Tr N U a (P, 0) is not always positive hence ( E ) (cr) can be ±1. Such fluctuations causes significant cancellations in 
the denominator of Eq. [TT] and renders the method ineffective due to large statistical uncertainties in (ft). In the 
literature, this problem is referred to as the Monte Carlo sign problem. For any Hamiltonian (Eq. [3l) with all V a < 0, 
h a are always time- reversal invariant, since all s a are real (Eq. [9J. As was shown by Lang et al. 15], time-reversal 
invariance of h a implies that the eigenvalues of the matrix come in complex-conjugate pairs which, in turn, ensures 
that the grand-canonical partition function TrU a is positive definite. In the canonical ensemble, projections on even 
number of particles always preserve the good sign as long as the grand canonical partition function is positive definite. 
However for systems with odd-number of particles, projections on to an odd number of particles usually reintroduces 
the sign problem at large values of f3 even when the grand canonical partition function is positive definite. 



IV. SIGN PROPERTIES OF THE ZERO-RANGE INTERACTION 



We now consider the zero-range interaction in the J-coupled representation which is discussed in full detail in 
ap pen dix [Al We write the two-body hamiltonian in the so-called pairing (or particle-particle) decomposition as 
[3; lla ] , which for the zero-range interaction gives 



ff 2 = \ E E V j( ab > cd)A\ M (ab)A JM (cd), (12) 

abed JM 
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where the antisymmetrized matrix element is given in eq. IA8I of appendix [XJ The pair operators are 

A jAi( ab ) = (3arn a j b m b \JM)a] bmb a] ama . (13) 

We now define the combined indices i = (ab) and j — (cd) to write Vj(i,j) which is a symmetric matrix. The goal 
is to diagonalize this latter matrix and look at the signs of the eigenvalues. As is demonstrated in [l5| , the interaction 
will have no sign problems when all eigenvalues of Vj(i,j) are negative. Obviously the problem splits into blocks of 
given J so we work in a fixed J subspace. 

The crucial observation is that Vj(a6, cd) can be factorized in the following way. First define the quantity 

fj(ab) = ^ { -i)W^[i a ] M j a ][ Jb] j£ i« J j ^ h e ie^M rRnJa{r)Rnbh{r)j (14) 

where e 2ie — sgn(Vo). Notice that fj(ab)e~ is a purely real number. In terms of the combined indices we now have 

Vj(i,j)= / drfj(i)fj(j). (15) 
Jo 

Since this matrix is real symmetric, there is a basis of orthonormal eigenvectors. Let us denote this basis u k and the 
corresponding eigenvalues A fc . The dimension is given by the number of pairs in the given model space that can couple 
to total angular momentum J. Consider now for a given k the product (u ) T Vju , where T denotes the transpose. 
Inserting the explicit form of Vj we have 

u T Vju = J2 i:j u k (i)Vj(i,j)u k ti) = 

I °°dr[Y,Jj(i)u k (i)] 2 = \\ (16) 

where the last equality follows from the eigenvalue equation and the fact that u k is normalized. We thus see that the 
eigenvalues are equal to some real number squared times a phase e 2l ° = sgn(Vo). Therefore, the sign of Vq is also the 
sign of the eigenvalues. We thus have the result that any attractive zero-range interaction (Vb < 0) will have no sign 
problem, whereas the repulsive (Vb > 0) case can never give a positive-definite path integral. 

The simple form of Vj(i,j) always us to prove some further properties of its spectrum. Define (for fixed J not 
shown) the row vector / = [/1/2 ■■■/«], where n counts the pairs, such as to fulfill Vj = J drf T f. Now pick a row 
vector orthogonal to / so that fg T = 0. Then we see that Vjg T — J drf T fg T — f drf T (fg T ) = 0, thus all vectors 
orthogonal to / are in the null-space of Vj. We therefore have only one non-zero eigenvalue for each J and n — 1 
eigenvectors with zero eigenvalue. The sole non-zero eigenvalue has the value / drff T and the eigenvector f T . We 
thus see that the zero-range pairing interaction has a very simple structure after diagonalization. 

As mentioned, the above proof was carried out in the so-called pairing decomposition with the operators A JM (ab) 
and AjM(cd). In many nuclear applications of the methods, the calculations are carried out in the density decom- 
position [7j. However, the exact path integral is independent of the particular representation and the above result 
will still hold. In particular, the change from pairing to density decomposition is in practice a re-coupling of the 
angular momenta involved (and a change of the one-body terms that we are not concerned with). Since re-couplings 
corresponds to changes of basis the result for the eigenvalues still holds. We have also done explicit numerical checks 
of this fact and confirmed the general statement. 

The good sign properties of the zero-range pairing rested on the fact that it could be factorized, which is more 
commonly referred to as separability of the zero-range interaction. A non-zero range interaction would not have this 
property and positive eigenvalues with associated sign problems can be expected. 

V. CONNECTION TO NUCLEAR PAIRING 

Many nuclear pairing studies consider only pairs of particles in time-reversed states with an attractive zero-range 
interaction of constant magnitude g < 0. This effectively means that the radial part of the matrix element in eq. IA6I 
is replaced by a constant. For the time-reversed states one has l a = h and l c = Id- A special case of the proof above 
thus shows that such a restriction preserves the good sign. 

The restriction to this interaction is inspired by the BCS approach to the pairing of electrons in conventional low 
temperature superconductors. The BCS theory was applied to nuclei only one year after its invention [16(. In the 
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general case, the gap will be state-dependent, Aj, where i denotes single-particle levels in the nuclear mean- field (as 
arising for instance in a Hartree-Fock calculation). However, the matrix element in the gap equations is the one in 
eq. lA6l with l a — 4 and l c — Id (see f.x. [13] )■ This approach is sometimes generalized into the Hartree-Fock-Bogoliubov 
(HFB) theory (also known as Bogoliubov-de Gennes in the atomic community), where pairing enters explicitly into 
the self-consistent equations for the optimal slater-determinant state [18j |. However, since the exact form of the pairing 
interaction in nuclei is not known, these studies are mostly phenomenological, and residual interaction derived from 
realistic nucleon-nucleon interactions cannot reproduce the magnitude of the experimentally observed gap [17[ . 

Notice that the SMMC described above would reduce to a HFB approach if the auxiliary field was determined 
self-consistently. However, as it is integrated out, the SMMC is in principle exact in the model space chosen and thus 
goes beyond HFB theory. It would be very interesting to see the differences between the full interaction and that 
acting only in time- reversed states within the SMMC approach. Such studies are currently under way @. 

VI. OUTLOOK 

Studies of small two-component fermi systems in tight traps are currently being pursued experimentally 0] and 
theoretically 0, Q . Here we demonstrates how the mapping of the atomic system to a nuclear physics problem can 
be done, and prove the very important feature that there is no sign problem with the typical choice of a zero-range 
interaction. As we have discussed, the atomic interaction between two internal hyperfine states is more general 
than the typical pairing force used in many nuclear investigations, and it was therefore not a priory clear that the 
transfer to nuclear shell-model Monte Carlo would be sign problem free. Of course one could also use large-scale shell- 
model diagonalization codes on the atomic problem. However, here the number of configurations grows exponentially 
with model space size, whereas in the SMMC it only grows quadratically Q. Of course it might be possible to do 
reasonable truncations and use effective interactions to reduce the size of the problem. Since low-dimensional systems 
are currently under intense study in atomic physics, this might also help alleviate to reduce the matrix problem in 
size. However, at the moment the shell-model Monte Carlo method seems to be the only tractable approach in the 
full three-dimensional case. 
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APPENDIX A: THE ZERO-RANGE FORCE IN J-COUPLED REPRESENTATION 

The matrix element used in the J-coupled representation for nuclear problems is given by 

V.,(ab,cd) = (kMrl) x 1> h ft)] JM \Vft,mi>i.(fi) * ^)] JM >, (Al) 

where a, 6, c, d have angular momentum j a ,jb,ja and jd, all of which are coupled from the external and internal 
quantum numbers. Notice that Vj is independent of the total projection M (this can be seen by applying the Wigner- 
Eckart theorem) . In a typical nuclear model there is some external mean-field (here assumed spherical) with quantum 
numbers (nlmi) and the internal spin half quantum numbers (|m s ). Coupling these through j = l + s gives j = Z± 1/2 
for the single-particle quantum numbers. 

As discussed in the text, in the atomic case with zero-range interactions the internal state must be a spin singlet; 
\S = 0, M s = 0). In order to ensure that the interaction in the J-scheme respects this, it is necessary to change the 
coupling scheme above from j a and jb coupled to J into l a and Zj coupled to L and s a and Sb coupled to S. This is 
done through the standard techniques of angular momentum [l9[ : 

\(l a s a )ja, (hs b )jb, JM) = Y,[L][S\\j a ]\j b ] \ h s b jl 1 \(l a l b )L, ( Sa S b )S, J M) . (A2) 

L,S [L S J) 
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Here we are interested in the s a = St, = 1/2 case, and, since the interaction contains a projection onto spin singlet 
states, only need the S = component of this transformation. Using a reduction on the 9j symbol the projection 
can be written 

Ps=o\m)ja, (hh)jb, JM) = E L {L}[ja}[jb] I i I % \ \(Ub)L, (||)0, JM) = 

[L J J 

( _ 1)i+ j a+ *+«.-i/ 2 I&^l i« jj * £J |(W 6 )L, (11)0, JM), (A3) 

where Ps=o = (1 — <T\ ■ <T2)/4 is the projection onto the spin singlet state. The remaining zero-range interaction of 
course only acts on the external quantum states, thus we have to evaluate matrix elements between coupled states 
with operators acting on only one of the degree of freedom. Since the spin part is trivial for singlets we simply have 
the result (keeping both L and J for clarity even though L = J) 

((l a l b )L, (||)0, JM\V(n - f 2 )|(Ui)£, (ii)0, JM) = 

((U^LMlVin - r 2 )\(l c l d )LM), (A4) 

where we have explicitly indicated the orbital angular momenta of all states involved. For the zero-range interaction 
V(f\ — r-z) — Vo5(f\ — F2), the latter matrix element can be found in many textbooks (see for instance (l9j ) and is 
given by 

(hhJMlVoSin-^W^J'M') = 
Sj,j>5 MlM ,[h][l2][im]$ o q) (o ^C^^ixWi^faWiZniJiCO^Cr). (A5) 
We can now insert all these formulae into eq. IA1I to get an expression for the J-scheme interaction: 

Vj(ab,cd) = ^ L (-i)'°+'c+2ja+2j e +j t +j d -i IMMMij] J« 
L Jc ^ ((/ a Z b )£0|U(ri - r 2 )\(l c l d )L0) = fj°lNWN J J Ja Jb\ J J 



1 / / f \\i.aib)^v\v V'l ~ '2;|V'ctd;^o/ - y-i) 2 \1 1 I [U / 7 

9 M *c J ^ 2 & a J I 2 a c 

(o o o) (o o) ^ir drr2R ^Ar)Rn blb (r)Rnj c (r)Rn d i d (r), (A6) 



where the second equality comes from using the formula in eq. IA5I Notice that the phase can be written with lb + Id 
instead of l a + l c since l a + l b + l c + U is even due to the restrictions from the Clebsch-Gordon coefficients. This is 
the general interaction in the spin singlet state and l a , lb, l c , and Id can in general be different as long as they couple 
pairwise to L — J. For the pairing interaction in the time-reversed states discussed in the text, we have l a = lb and 
l c = Id- 

The formula above explicitly shows that l a + lb + J and l c + + J must be even. However, the multipole expansion 
used to arrive at this expression implicitly requires that also l a + l c + J and lb + Id + J be even. Notice also that 
the factor l a + l c means that pairing across two opposite parity major shells can be repulsive for Vq < 0. This is 
well-known in nuclear pairing studies Q. 

As discussed in Q, the physical matrix elements used in the nuclear shell model must be antisymmetrized. This 
can be achieved by using the definition 

Vf(ab, cd) = -. 1 

[Vj{ab, cd) - (-iy*+^- J Vj(ab, dc)] . (A7) 
However, as one can easily show, Vj(ab,dc) = —(—l) : >' :+: > d ~ J Vj(ab,cd). Thus we have the simple result 

Vj A (ab, cd) = 2 =Vj(ab, cd). (A8) 

V(l + o ab ){l + d cd ) 
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This is not surprising since we argued that only the antisymmetric S — spin-singlet component should have non-zero 
matrix elements. We have effectively enforced the Pauli principle in this manner. 
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